Accepted by Astrophysical Journal Letters, 26 June 2002 

Angular momentum extraction by gravity waves in the Sun 

Suzanne Talon 

Departement de Physique, Universite de Montreal, Montreal PQ H3C 3J7 
CERCA, 5160, boul. Decarie, suite 400, Montreal PQ H3X 2H9 

Pawan Kumar 

Dept. of Astronomy, University of Texas, Austin, TX 78713 
Jean-Paul Zahn 
LUTH, Observatoire de Paris, 92195 Meudon, France 

Abstract 

We review the behavior of the oscillating shear layer produced by gravity waves be- 
low the surface convection zone of the Sun. We show that, under asymmetric filtering 
produced by this layer, gravity waves of low spherical order, which are stochastically 
excited at the base of the convection zone of late type stars, can extract angular mo- 
mentum from their radiative interior. The time-scale for this momentum extraction in 
a Sun-like star is of the order of 10 7 years. The process is particularly efficient in the 
central region, and it could produce there a slowly rotating core. 

Subject headings: hydrodynamics — stars: late-type — stars: rotation — Sun: interior 
— waves 

1. Introduction 

Angular momentum transport by gravity waves has received much attention recently. While it 
has first been considered as a key mechanism in the tidal interaction of binary systems (Zahn 1975; 
Goldreich & Nicholson 1989), it has since been invoked also as an efficient process of momentum 
redistribution in single stars (Schatzman 1993; Kumar & Quataert 1997; Zahn, Talon & Matias 
1997). The mechanism proposed was similar to that acting in binary stars, with synchronization 
proceeding gradually inwards. However, as pointed out by Gough & Mclntyre (1998) and Ringot 
(1998), the treatment of angular momentum extraction as presented by these first studies was 
incorrect, and it became clear that the actual properties of wave transport were far more complex 
than originally anticipated. 
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Fig. 1. — Angular momentum luminosity integrated over 0.1/uHz as a function of order I and 
frequency. 

Wave properties have been examined further by Barnes, McGregor & Charbonneau (1998) who 
showed that magnetic fields transform pure gravity waves into gravito-Alfven waves, thus modifying 
their damping behavior. They showed that strong magnetic field may prevent certain waves from 
propagating. 

Kumar, Talon & Zahn (1999, hereafter KTZ) then demonstrated how gravity waves in inter- 
action with shear turbulence may lead to the formation of an oscillating shear layer just below 
the surface convection zone, analogous to the well studied quasi-biennial oscillation of atmospheric 
sciences (see e.g. Shepherd 2000 for details on the QBO and other wave-driven oscillations). 

Kim & McGregor (2001) further studied that oscillation in a simplified two waves model. They 
showed that, with only a prograde and a retrograde wave, and with the velocity fixed both at the 
top and at the bottom of the solar tachocline, the time-dependent behavior of the rotation profile 
goes from being periodic, to quasi-periodic and to chaotic as the viscosity in the shear region is 
slowly decreased. Their goal was to seek an explanation of a 1.3 yr variation, which Howe et al. 
(2000) claim to have detected in the tachocline. 

Here we want to reexamine the effect of waves in the deep solar interior. Indeed, while the 
high degree waves are damped very close to the base of the convection zone, thus leading to rapid 
oscillations that could explain rotational velocity variations in the tachocline (Kim & McGregor 
2001) or be related to the solar cycle (KTZ), the low degree waves that are also excited stochastically 
by the convective motions, although with a somewhat lower efficiency (see Kumar & Quataert 1997 
and Zahn, Talon & Matias 1997) will have, over a much longer period, an effect on the deep interior. 

In this letter, we wish to demonstrate how these waves, in conjunction with shear turbulence, 
indeed achieve momentum extraction in the deep interior. In §2, we explain heuristically how that 
process takes places and in §3, we show the results of a long numerical simulation applied to the 
solar case. 



3 



an (,uHz) 

0.02 


-0.02 

0.02 


-0.02 

0.02 


-0.02 

0.7 0.705 0.71 0.7 0.705 0.71 0.7 0.705 0.71 

Fig. 2. — Oscillating shear layer below the surface convection zone. The dotted line shows the initial 
rotation profile. With the surface rotating slower than the core, a prograde shear layer is initially 
formed, followed by a retrograde one. When the shear becomes too intense, turbulent viscosity acts 
to merge the prograde layer with the convection zone, leaving behind the retrograde layer. A new 
prograde layer forms behind, and migrates towards the convection zone when the retrograde layer 
is absorbed. The cycle then resumes. 

2. Shear layer and momentum extraction 

The formation of a shear layer is easy to describe, using a simple two waves model. When both 
a prograde and a retrograde wave travel in an accelerating region (i.e. where Q increases inwards 
with depth), the prograde wave frequency diminishes, enhancing its damping; the opposite occurs 
for the retrograde wave (see KTZ for details). This creates a double peaked shear layer. The 
gradient in that layer rises until it becomes steep enough so that turbulence will start playing a role 
and will effectively merge the shear layer with the adjacent convection zone. The retrograde layer 1 
will then slowly rise toward the convection zone, and a prograde layer will form further inside the 
star. The retrograde layer is eventually also absorbed by the convection zone, and the cycle goes 
on. 

The oscillating shear layer will play an important role in filtering the (low degree) waves that 
will penetrate, and thus deposit their momentum in the deep interior. Indeed, when the prograde 
peek is larger, it is the retrograde waves that will penetrate and deposit their negative momentum in 
the interior, and vice versa if the retrograde peek dominates. If there is no initial differential rotation 
in the interior, the shear layer will oscillate symmetrically, and on average no net momentum flux 
will occur. However, if the surface convection zone is slowed down, as in the case of the Sun, there 
exists a bias in the formation of the shear layer, which will favor the prograde peek (the interior 




lr The retrograde layer is the one that rotates more slowly than the convection zone, while the prograde layer rotates 
more rapidly. 



- 4 - 




Fig. 3. — Global evolution of the rotation profile over large time periods. Efficient momentum 
extraction is visible in a characteristic time-scale of 10 5 years. 

then rotating faster than the surface), thus letting a larger amount of negative momentum being 
transported to the deep interior. The positive momentum deposited in excess in the shear layer will 
then be transfered back in the surface convection zone by turbulence. The net effect will thus be to 
extract momentum from the radiation zone. 



3. A numerical simulation 



While it is clear that the filtering bias introduced by the slope of the rotation profile below 
the convection zone leads to momentum extraction, we further need to know where exactly this 
extraction takes place, and how efficient it is. To answer this question, we performed a numerical 
simulation covering a very long time period. There are two key ingredients in such a simulation. 
Firstly, one needs to include a whole variety of radial orders to describe both the waves responsible 
for the formation and oscillation of the shear layer and the deeply penetrating waves. Furthermore, 
retaining only two wave orders has the effect of depositing momentum too locally, compared to 
having a more complete spectrum. We thus included spherical orders £ = 1, 5, 9, 57, to capture 
a wide variety of damping behaviors. We also considered a small frequency range u> = 0.6 to 
1.0 fiHz. Higher frequency waves will have little impact since the power spectrum we are using 
in these simulations falls off as u>~ 4 ' 5 (cf. Eq. 5 and KTZ), while lower frequency waves will be 
damped very close to the convection zone (cf. Eq. 3), not undergoing differential damping. The 
second key ingredient is to follow both short time-scales (~ 10 years), characteristic of the shear 
layer oscillation, and long time-scales (~ 10 7 years), characteristic of the spin-down rate. This 
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Table 1: Damping length d w (scaled by the solar radius) at 0.6 /jHz 



t 


m 


d a 


d b 


1 


1 





616 


0.461 




-1 






0.661 


5 


5 





084 


0.036 




-5 






0.499 


9 


9 





033 


0.017 




-9 






0.378 


13 


13 





018 


0.010 




-13 






0.256 


17 


17 





012 


0.007 




-17 






0.061 


21 


21 





008 


0.005 




-21 






0.024 



a In solid body rotation 

6 With the initial rotation profile (c/. Figures 2 & 3) 

combinaison requires the use of a conservative numerical scheme. 

We thus calculated the evolution of angular momentum in the radiative interior by the combined 
effect of turbulence and gravity waves: 



^ dt r 2 9r 



pv t r 

or 



where p is the density, fl the angular velocity and u t the turbulent viscosity (cf. Talon 1997). The 
local momentum luminosity Cj{r) is integrated over the whole wave spectrum, each wave being 
damped by radiative processes and turbulent viscosity (see Eq. 7) 

= Cji > m ( rzc ) exp t~ r ( r ' a ' ^ ( 2 ) 

cr,£,m 

and where the local damping rate takes into account the mean molecular weight stratification 

? r Tc NN 2 / N 2 \ 5 r\r 



T 



where N 2 = Nj, + N 2 is the total Brunt-Vaisala frequency, Nj, is its thermal part and N 2 is due to 
the mean molecular weight stratification (cf. Zahn et al. 1997). a(r, m) = to — m [£l{r) — fl zc ] is the 
local (Doppler shifted) frequency, and u is the frequency in the reference frame of the convection 
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zone. The damping length d w defined by T(d UJ ) = 1 depends strongly on frequency u and, in 
presence of differential rotation, on the signed value of the azimuthal wavenumber m, as shown in 
Table 1. The evolution equation (1) is solved using a decoupled scheme, in which momentum is 
first deposited by waves and diffusion is then treated using a perfectly conservative finite element 
method. After 10 7 time steps of 1 year, the error on global momentum conservation is ~ 10 -5 . 
This equation is solved throughout the radiation zone. The upper boundary condition expresses 
the conservation of momentum of the star as a whole, with the convection zone rotating as a solid 
body. On their way to the core, the waves are reflected when their local frequency a equals the 
Brunt-Vaisala frequency 

a(r,m) = N. (4) 

They then deposit more momentum as they travel back to the convection zone. Some waves do 
travel all the way to the core. Since that region contains very little momentum, it is easily spun 
down (see Figure 3). Waves reaching the center most region (below r = 0.01 R & ) are reflected back; 
this innermost core has been taken to rotate as a solid body. 

Wave generation is due to Reynolds stresses at the base of the convection zone, and to model 
its amplitude we followed here the description of Goldreich et al. (1994) (see also KTZ). The energy 
flux per unit frequency Te is then 

d &\ 2 ,„,„,-,, f d£ h 



v 3 r4 

xexp[-^ + l)/2r 2 ] i + ( ^ )15/2 , (5) 

where £ r and [£(£ + l)] 1 / 2 ^ are the radial and horizontal displacement wavefunctions which are 
normalized to unit energy flux just below the convection zone, v is the convective velocity, L is 
the radial size of an energy bearing turbulent eddy, tl ~ L/v is the characteristic convective time, 
and hu is the radial size of the largest eddy at r with characteristic frequency of co or greater 
(h u = Lmin{l, (2u;tl) _3 / 2 }). The gravity waves are evanescent in the convection zone, the region 
where they are excited. The above equation was derived under the assumption that the turbulence 
spectrum is Kolmogorov and it ignores wave excitation resulting from convective overshooting. The 
momentum flux per unit frequency Tj is then 

777 

Tj{m,l,u) = -T E {l,w). (6) 

id 

In the case of the Sun, there is a whole spectrum of waves, in frequency ui and spherical 
order £. Their range is given by to < N c , where N c is the Brunt-Vaisala frequency at the base 
of the convection zone 2 , and by 1 < £ < £ c where £ c ~ 60 is the spherical order characterizing 



2 This frequency depends on the amount of overshooting that is considered. However, since the power spectrum 
scales as u>~ 4 ' 5 (cf. KTZ), its exact value is not crucial. 
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convection (it corresponds to one pressure scale-height). Assuming that all azimuthal orders m are 
equally excited, we find that most of the momentum is carried by the low frequency waves. There 
is significant momentum luminosity in the low order waves that do penetrate deep in the interior 
(cf. Figure 1 and Table 1). 

For this first exploration, turbulent viscosity has been taken proportional to the radiative 
viscosity (f ra d), admittedly a very crude prescription: 

^ = 10 5 z, rad = 10 5 ^ (7) 

where x 1S the thermal conductivity and c the speed of light. The total energy flux in the considered 
waves is 6.3 x 10 29 erg/s, and the total integrated momentum luminosity is 1.7 x 10 36 gcm 2 /s 2 . This 
turbulent viscosity operates both on the mean rotation profile and on the waves (cf. Eq. 3 and 1). 

This simulation is performed in one dimension and we only consider the effect of spherically 
integrated waves. We thus neglect the latitudinal differential rotation just below the convection 
zone, as well as the influence of the Coriolis force on wave propagation and wave excitation. 

The results are displayed in Figures 2 and 3. The first shows the shear layer oscillation cycle, 
which occurs here with a period of ~ 300 years. The second shows the evolution of the rotation 
profile in the deep interior over a much longer period. Momentum extraction is related to the 
asymmetric filtering by the double shear layer, which is determined by the amount of differential 
rotation initially present below the convection zone. In our simulation we let O increase by 0.15 fiHz 
down to r = O.6Oi? . This profile, applied to the wave spectrum we consider here, yields a differential 
momentum flux below the shear layer of 2x 10 33 gcm 2 /s 2 . In 5x 10 4 years, this leads to an extraction 
of 2.5 x 10 45 gcm 2 /s, in good agreement with the numerical simulation which shows an extraction 
of 4 x 10 45 gcm 2 /s. 

As the simulation progresses, the radiation zone decelerates and the convection zone acceler- 
ates 3 ; this leads to less differential filtering and therefore momentum extraction gradually slows 
down. The simulation shows that waves are very efficient in extracting momentum from the deep 
core, in a time short compared to the main-sequence life time. At later times, the net wave flux 
decreases, allowing turbulence to smooth the rotation profile. But the turbulent transport will also 
taper off, due to the build-up of a stabilizing helium gradient. Therefore it is plausible that in the 
present Sun the central region is rotating slower than the rest of the radiation zone, as suggested 
by certain helioseismic inversions (Elsworth et al. 1995; Corbard et al. 1997). If this trend is 
confirmed, it would imply that wave transport is more efficient than magnetic torquing, which has 
also been invoked for the extraction of angular momentum (cf. Mestel & Weiss 1987, Charbonneau 
& MacGregor 1993). 



3 In an evolving solar model, the surface convection zone would constantly be spinned down by a magnetic torque 
and would thus probably never accelerate. 



- 8 - 



4. Discussion 

In this letter, we have shown how gravity waves can extract angular momentum from the 
radiative core of a solar-like star through differential wave filtering, sufficiently to establish an almost 
uniform rotation profile in about 10 7 years. In our model, the presence of turbulent transport is 
essential to yield the proper result; indeed, it is required both to produce the QBO like oscillation, 
and to smooth the excessive differential rotation produced by waves in the central region. 

The evolution of latitudinally averaged angular momentum in the Sun would then proceed as 
follows. The magnetic torque applied to the convection zone tends to slow it down, establishing 
a negative gradient of angular velocity, which induces the bias between prograde and retrograde 
waves. Due to the combined action of wave transport and of turbulent viscosity, the velocity 
gradient tends to diminish. In our simulation it disappears altogether, but in the real Sun, it would 
be maintained at some level due to the competition between the extraction of angular momentum 
from the radiative interior, and its loss by the solar wind: the faster the spin-down, the steeper the 
gradient, the stronger the extraction. If the extraction of angular momentum by gravity waves is as 
efficient as we estimate, it will adjust itself such as to provide just the flux of angular momentum 
which is lost by the wind, with the interior profile of angular velocity being rather flat, except near 
the center and at the top of the radiation zone. 

The presence of a large toroidal magnetic at the base of the convection zone could somewhat 
lengthen the time-scales mentioned here. Indeed, the magnetic field required to prevent wave 
propagation at 0.6 /uHz is about (3 x 10 5 /£) G (c/. KTZ). If fields of a strength as high as 10 5 G are 
present in that region (c/. Fan et al. 1993, Caligari et al. 1995), high i modes could be prevented 
from propagating, modifying the evolution of the double shear layer. However, since such large 
fields occur in localized areas, the global effect would rather be one of decreasing the wave flux and 
thus, simply increase the time scale estimated here. 

This initial work must be pursued in many ways. Firstly, one would like to assess the efficiency of 
this mechanism in an evolving solar model, with magnetic spin-down, in which meridional circulation 
would also be taken into account in order to verify the scenario we proposed earlier. 

Other improvements should be made on the physical description. The effect of the Coriolis 
force must be included, which could require to step up to two-dimensional simulations. Finally, a 
more realistic prescription for the turbulent transport must be implemented, before we can conclude 
that the solar core is rotating significantly slower than the rest of the radiation zone. 
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